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Abstract 

Background: Small non-coding RNAs (ncRNAs) are innportant regulators of gene expression in eukaryotes. 
Previously, only microRNAs (nniRNAs) and piRNAs have been identified in the silkwornn, Bombyx mori. Furthermore, 
only ncRNAs (50-500nt) of intermediate size have been systematically identified in the silkworm. 

Results: Here, we performed a systematic identification and analysis of small RNAs (18-50nt) associated with the 
Bombyx mori argonaute2 (BmAgo2) protein. Using RIP-seq, we identified various types of small ncRNAs associated 
with BmAG02. These ncRNAs showed a multimodal length distribution, with three peaks at ~20nt, ~27nt and 
~33nt, which included tRNA-, transposable element (TE)-, rRNA-, snoRNA- and snRNA-derived small RNAs as well as 
miRNAs and piRNAs. The tRNA-derived fragments (tRFs) were found at an extremely high abundance and 
accounted for 69.90% of the BmAgo2-associated small RNAs. Northern blotting confirmed that many tRFs were 
expressed or up-regulated only in the BmNPV-infected cells, implying that the tRFs play a prominent role by 
binding to BmAgo2 during BmNPV infection. Additional evidence suggested that there are potential cleavage sites 
on the D, anti-codon and TipC loops of the tRNAs. TE-derived small RNAs and piRNAs also accounted for a 
significant proportion of the BmAgo2-associated small RNAs, suggesting that BmAgo2 could be involved in the 
maintenance of genome stability by suppressing the activities of transposons guided by these small RNAs. Finally, 
Northern blotting was also used to confirm the Bombyx 5.8 s rRNA-derived small RNAs, demonstrating that various 
novel small RNAs exist in the silkworm. 

Conclusions: Using an RIP-seq method in combination with Northern blotting, we identified various types of small 
RNAs associated with the BmAgo2 protein, including tRNA-, TE-, rRNA-, snoRNA- and snRNA-derived small RNAs as 
well as miRNAs and piRNAs. Our findings provide new clues for future functional studies of the role of small RNAs 
in insect development and evolution. 
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Background (piRNAs). Mature miRNAs are a class of ~22nt small 

Small non-coding RNAs (ncRNAs) are important regu- RNA molecules that bind to complementary sequences 

latory molecules that regulate the expression of target in target genes and regulate gene expression by modu- 

genes by affecting mRNA stability, protein translation lating mRNA translation or stability [1]. miRNAs have 

and chromosome modification. There are three main been implicated in many cellular processes, including 

classes of small RNAs: microRNAs (miRNAs), small- immune mechanisms of biological molecules, stem cell 

interfering RNAs (siRNAs) and Piwi-interacting RNAs differentiation, cell signal transduction, tumorigenesis, 

neuronal development and apoptosis [1,2], siRNAs are a 
class of double-stranded RNA molecules between 20 and 
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binding sites [3]. siRNAs are also known to play an 
important role in antiviral mechanisms and in shaping the 
chromatin structure of the genome [4-6]. Recently, many 
endogenous siRNAs were identified that function in tran- 
scriptional silencing [7-9]. piRNAs are a class of small 
RNAs associated with PIWI family proteins. The biogen- 
esis of piRNAs is related to a process called the Ping-pong 
cycle [10,11]. These small RNAs have been found mainly 
in stem cells and germ cells, and they are thought to play 
an essential role in germline development, stem cell 
renewal, transposon silencing and epigenetic regulation 
[12-15]. Due to the recent development and application of 
high-throughput sequencing, many novel small RNAs 
have been identified, including tRNA-, rRNA-, snoRNA- 
and snRNA-derived small RNAs [16-21]. tRNAs are not 
very stable and, under stress, can be sequence-specifically 
cleaved into different sized fragments termed "tRNA- 
derived fragments," "tRFs" or "sitRNAs." These fragments 
could be a novel class of small RNAs. Because they associ- 
ate with the Ago2 protein, it seems likely that they play a 
role similar to siRNAs in RNAi silencing [16,17]. 

The silkworm, Bombyx morU is an important economic 
insect and has become a model organism for the study 
of the biochemistry, molecular genetics and genomics 
oi Lepidoptera insects [22,23]. Previous studies on small 
ncRNAs in the silkworm have focused on miRNAs and 
piRNAs. Our group was the first to provide a large-scale 
identification of miRNA genes in Bombyx mori [24]. 
Recently, a large number of new miRNAs have been 
identified using deep sequencing [25-29] and their 
potential functions were predicted through an analysis 
of expression profiling performed at various developmental 
stages and in various tissues of the silkworm [26,30-33]. 
The miRBase database (http:www.mirbase.org/. Release 
19.0) has published 567 Bombyx miRNAs. piRNAs have 
also been well characterized in the silkworm. Kawaoka 
et al analyzed the biogenesis of piRNAs, which could 
exert an important genomic defense against transposons 
in the silkworm genome [34-40]. However, less work on 
siRNAs in the silkworm has been performed, and only 788 
potential transposable element (TE) -associated siRNAs 
have been identified by deep sequencing techniques [18]. 
In addition, intermediated- sized ncRNAs (50-500nt) have 
been systematic identified in the silkworm, including 141 
snoRNAs, six snRNAs and 38 unclassified ncRNAs [41]. 
Based on the recent identification of an increasing number 
of small RNAs, it seems likely that many novel small 
RNAs remain to be discovered in Bombyx mori. 

To execute their functions, small RNAs must be incor- 
porated into the RNA-induced silencing complex (RISC). 
Once incorporated, the small RNAs guide RISC to its 
complementary targets and thereby regulate the expres- 
sion of the targets at the transcriptional or translational 
levels [42,43]. The Argonaute (Ago) protein is a central 



protein component of RISC and exerts specialized func- 
tions on the miRNA and RNAi pathway [44-46]. The 
mammalian Argonaute family has four members, but only 
Ago2 possesses slicer endonuclease activity [46,47]. The 
isolation of RNAs associated with Ago proteins has be- 
come an important means for the discovery and functional 
assessment of novel small RNAs [18,21,48-52] and for the 
target genes of miRNAs [53,54]. Bombyx mori Argonaute2 
(BmAgo2, GenBank accession: NM_001043530.2) belongs 
to the Ago family and is an ortholog of Drosophila 
Argonaute2, which contains the conserved amino acid 
residues D965, D1037 and HI 173. These conserved 
residues are critical for the nuclease activity of Ago2. 
Previous reports have shown that in silkworm infected 
with Bombyx mori nucleopolyhedrovirus (BmNPV), 
BmAgo2 expression is up-regulated, which could be 
related to the RNA silencing machinery involved in 
DNA virus infection in insects [55,56]; however, this 
mechanism will require further study. 

Ago proteins are key components of the siRNA and 
miRNA pathway and are indispensable binding proteins 
for the function of many other small RNAs. Therefore, the 
isolation of Ago-associated small RNAs is an important 
approach for identifying functional small RNAs [18,19,21]. 
In this study, we extracted the total small RNAs (18-50nt) 
that associated with BmAgo2 protein using the RNA 
immunoprecipitation (RIP) method. Subsequent deep 
sequencing, bioinformatics analysis and Northern blot- 
ting were used to identify various types of small RNAs 
associated with the BmAgo2 protein, including tRNA-, 
TE-, rRNA-, snoRNA- and snRNA-derived small RNAs as 
well as miRNAs and piRNAs. Further analysis revealed 
that these small RNAs possess novel characteristics. 

Results 

RIP of BmAgo2 from BmN cells Infected with recombinant 
BmNPV virus 

Small RNAs and their targets bind the Ago-containing 
RISC complexes, in which the Ago proteins form stable 
Ago ribonucleoproteins that can be biochemically ana- 
lyzed [53,57,58]. The Ago-protein-binding small RNAs 
can be isolated by RIP [59,60]. In a previous work, 
BmAgo2 was fused with a HIS tag and was successfully 
expressed using the Baculovirus Bacmid system harbor- 
ing the iel promoter enhanced with a hr5 enhancer 
[61]. The recombinant viruses were then harvested at 
20 hrs post infection, and HIS-BmAgo2 could be 
detected at a high level by Western blotting with a HIS 
monoclonal antibody (Additional file 1: Figure SI). The 
HIS monoclonal antibody (mouse anti-(his)6, Roche) was 
used to immunoprecipitate HIS-BmAgo2-containing RISC 
from the total cell lysate of the infected BmN cells. The 
approximately 120 kDa HIS-BmAgo2 was identified by 
Western blotting in the total cell lysate and HIS-BmAgo2 
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IP fraction but was absent in the IP fraction of the nega- 
tive control (Figure lA). The co-immunoprecipitated 
BmAgo2-bound RNAs were extracted and analyzed by 
PAGE. Interestingly, the RNA collected via the HIS- 
BmAgo2-specific monoclonal antibody pull-down showed 
a much more dense RNA smear than the total RNAs of 
the BmN cells, ranging from 18 bp to 50 bp (Figure IB). 
In this region, three visible bands of ~20nt, ~27nt and 
~33nt accumulated, indicating that small RNAs may be 
recruited by the BmAgo2 complexes and may function 
together with the BmAgo2 protein. The negative control 
IgGl from the mouse did not immunoprecipitate visible 
amounts of RNAs (Figure IB), indicating that the 
BmAgo2-associated RNA signal was specific. Taken to- 
gether, our results demonstrate that BmAgo2 can bind 
many small RNAs that can be biochemically isolated. 

Sequencing and analysis of a BmAgo2-associated small 
RNA library 

Ago proteins have an important role in small RNA 
pathways and mediate the interaction between small 
RNAs and their targets. The resolution of Ago- 
associated small RNAs showed a significant landscape 
of Ago proteins and their binding to small RNAs 
[18,19,21,48-50]. To characterize the small RNAs that 
associate with the BmAgo2 protein in Bombyx mori, 
the small RNA population associated with this protein in 
BmN cells was extracted from the Ago immunoprecipitated 



complex. Small RNAs between 18nt and 50nt were 
separated by PAGE and were then subjected to library 
construction and deep sequencing. 

The high throughput sequencing yielded a total of 
813,702 unique reads, representing 11,691,441 high- 
quality total reads of between 1 and 5,731,905 copies. The 
reads were then mapped to the Bombyx mori genome. Of 
the 813,702 unique reads, only 538,567 mapped to the 
reference genome when mismatches were not allowed, ac- 
counting for 66.19% of the total unique reads and showing 
no preference for the distribution between the chromo- 
somes; however, when a maximum of two mismatches 
were allowed, 722,610 unique reads mapped to the refer- 
ence genome, accounting for 88.81% of the total unique 
reads, suggesting that many mutations, base edits or SNPs 
exist in these small RNAs (Additional file 1: Figure S2A). 
Thus, 11.19% of the total unique reads did not match the 
genome, which may also be due to possible gaps in the se- 
quenced Bombyx mori genome [62,63]. In addition, the 
reference genome was derived from male larvae and there- 
fore did not contain sequences from the W chromosome 
[62,63]. The sillcworm W chromosome was a source of 
small RNAs in Bombyx, such as the female-enriched 
piRNAs [36]. Of the 91,092 unmapped unique reads, 
9,206 mapped to contig sequences from the W chromo- 
some, which form a small portion of the whole W 
chromosome. Additionally, 218 unique reads derived from 
the BmNPV genome were identified. 
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Figure 1 RIP for the Isolation of BmAgo2-associated RNAs from BmN cells. (A) Co-lmmunoprecipitation of HIS-BmAgo2 using HIS 
monoclonal antibody. IP was validated by Western blotting using a HIS monoclonal antibody. BmAgo2 was successfully expressed in BmN cells 
infected with recombinant BmNPV virus (lane 1) and was pulled down when the HIS antibody was used (lane 3). Mouse IgGl was used as a 
negative control and indeed revealed no IP of BmAgo2 (lane 2). (B) Co-immunoprecipitated RNAs analyzed by PAGE. Compared to the total RNAs 
of the BmN cells (lane 3), the HIS monoclonal antibody pull-down RNAs showed an obvious RNA smear ranging from 20 bp to 50 bp (lane 1). The 
negative control mouse IgGl did not immunoprecipitate visible RNA amounts (lane 2). Lane M shows a DNA size marker. Asterisks indicate the visible 
bands for the ~20nt, ~27nt and ~33nt small RNA species, respectively. 
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We investigated the length distribution of the 
obtained BmAgo2-associated small RNAs. The reads 
showed a multimodal length distribution, with three peaks 
at ~20nt, --27^, and --SSnt (Additional file 1: Figure S2B). 
These sizes were consistent with the sizes of the bands 
separated by PAGE (Figure IB). The length distribution of 
the BmAgo2-associated small RNAs was different than 
that of the SjAgo2-associated small RNAs, which showed 
only one peak at 21nt [18]. The number of 32nt and 
33nt small RNAs in the BmAgo2-associated total small 
RNA population decreased significantly when the unique 
sequences were analyzed, but the number of ~20nt 
and ~27nt small RNAs did not change much. This 
finding suggests that the BmAgo2-associated ~20nt 
and ~27nt small RNAs comprised a large number of 
diversified sequences with low expression abundance, 
whereas the expression of the ~33nt small RNAs oc- 
curred at an extremely high level (Additional file 1: Figure 
S2B and C, Additional file 2: Table SI). Additionally, 
we observed a strong preference for the first 5 ' nucleo- 
tide to be "U" in the unique BmAgo2-associated small 
RNA sequences (accounting for 54.7% of the unique 
small RNA sequences) and "A" was similarly enriched 
at position 10 (Figure 2A), consistent with previous 
reports [19,34,35,64,65]. 

The obtained small RNAs were then classified and 
annotated using the bioinformatics pipeline described 
in the Materials and Methods. The results showed that 
"housekeeping" non-coding RNAs, such as tRNAs, 
rRNAs, snRNAs and snoRNAs, could produce BmAgo2- 
associated small RNAs. Thus, various types of small 
ncRNAs were identified among the BmAgo2-associated 
small RNAs, including TE-, tRNA-, rRNA-, snoRNA and 
snRNA-derived small RNAs as well as miRNAs and 
piRNAs (Figure 2B). In addition, the 218 BmNPV-derived 
unique reads might be novel BmNPV-encoded small RNA 
species. The unique small RNA species were dominated 
by TE-derived small RNAs, which accounted for 47.78% 
of the unique small RNAs, but the total small RNAs were 
overwhelmingly dominated by tRFs, which accounted for 
2.21% of the unique small RNAs and 69.90% of the total 
small RNAs. We also identified many small RNAs without 
clear structural features or functions, which we defined as 
unannotated small RNAs. The unannotated small RNAs 
comprised a large number of diverse sequences with low 
expression and accounted for 38.91% of the unique small 
RNAs and 9.00% of the total small RNAs. 

tRFs 

With the development of molecular cloning for the 
identification of cellular RNA fragments and deep se- 
quencing, a number of classes of small RNAs have been 
identified recently, including a novel class of small 
RNAs in eukaryotes. This group of molecules, termed 



tRFs, comprises RNA fragments derived from either 
tRNAs or tRNA precursors [16,17,20,66]. tRFs were the 
most abundant small RNA identified from the deep 
sequencing of the BmAgo2-associated small RNAs, 
accounting for 69.90% of the total small RNAs. How- 
ever, tRFs accounted for only 2.21% of the unique small 
RNAs, implying that an extremely high abundance but a 
low species number of tRFs are recruited by the 
BmAgo2 protein (Figure 2B). 

To further confirm the expression of the tRFs and to de- 
termine their size, we used Northern blotting to analyze 
the selected tRFs. The tRFs with read copy numbers of no 
less than 150 and redundancy removal were selected. In 
addition, the other selected tRFs were derived from the 
same host tRNA that generate the previously selected 
tRFs. Our results showed that Bombyx tRNAs could 
generate small RNAs, which were derived mainly from 
the 5' or 3' termini of the tRNA (termed the 5'tRF or 
3'tRF, respectively). There were two types of 5'tRF 
molecules. One was a ~33nt 5'tRF, which was located at 
the 5' termini of the host, tRNA and produced via 
cleavage of the anti-codon loop by a certain nuclease 
(Figure 3A and Additional file 2: Table S2). The other 
was a truncated 18nt 5'tRF molecule, which was pro- 
duced by further cleavage of the ~33nt 5'tRF at the D 
loop; this 5'tRF molecule could not be detected using 
Northern blotting but could be sequenced at a high read 
copy number (Bm-SerAGA-5 'T and Bm-LeuCAG-5 'T, 
Additional file 2: Table S2). The 3 'tRFs also comprised 
two types of molecules. One was an ~40nt small RNA 
that was located at the 3 ' termini of the host tRNA and 
was produced through the cleavage of the anti-codon 
loop (Figure 3B and Additional file 2: Table S2). The 
second was a 21nt 3 'trailer tRF small RNA located 
between the anti-codon and Tv|/C loop of the host tRNA 
(tRNA 3 'trailer) that was produced by further cleaving 
the --40^ 3'tRF at the Tv|/C loop. The 3 'trailer tRF 
could not be detected using Northern blotting but could 
be sequenced with a very high read number (Bm- 
MetCAT-3'T, Bm-GlnCTG-3'T, Bm-GlyCCC-3'T and 
Bm-IleTAT-3'T, Additional file 2: Table S2). A GTTC 
motif appeared in the 3 'termini of the 21nt 3 'trailer tRFs 
due to a conserved cleavage site at the end of the GTTC 
motif in the Tv|/C loop (Additional file 1: Figure S3). The 
asymmetric source of tRFs from the tRNAs showed that 
the tRFs are not random by-products of tRNA biogenesis 
and degradation but are instead an abundant and novel 
class of small RNAs. 

Interestingly, Bombyx tRFs seem to play a prominent 
role in the response of Bombyx mori to BmNPV. Many 
tRFs showed very high expression levels in the 
BmNPV-infected BmN cells (Figure 3). Among these 
tRFs, some exhibited obviously up-regulated expres- 
sion in the BmNPV-infected BmN cells, while others 
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Figure 2 Characterization of BmAgo2-associated small RNAs. (A) Nucleotide compositions of unique small RNAs. BmAgo2-associated small RNAs 
displayed a strong bias towards "U" at the 5' end (54.47%). (B) Classification and annotation of small RNAs using the bioinformatics pipeline described 
in the Materials and Methods. tRFs were recruited by the BmAgo2 protein at an extremely high abundance, and the unannotated small RNAs 
comprised many diverse sequences with a low expression abundance. Of the 1 10,157 piRNA candidates, 54033 were derived from TEs, accounting for 
6.64% of the total unique small RNAs. Thus, the actual percent of TE-derived sRNAs was 47.78% (41 .14% + 6.64%). sRNAs: small RNAs. 



were expressed only in the BmNPV-infected BmN cells during baculovirus infection. Importantly, almost all of the 
(Figure 3A). These results suggest that the tRFs are im- detected 5 'tRFs were expressed in only the infected BmN 
portant stress factors and might have an important function cells; thus these tRFs are most lil<ely produced under the 
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Figure 3 Identification of tRFs by Northern blotting. (A) 5'tRFs were identified by Nortliern blotting witli a lengtli of -33 bp. Tine majority of 
tine identified 5'tRFs were generated only in BmNPV-infected cells, implying a BmNPV-related role for these small RNAs. (B) 3'tRFs were identified 
by Northern blotting with a length of -40 bp. However, the majority of the identified 3'tRFs were generated in both the normal cells and in the 
BmNPV-infected cells. The host tRNAs were also identified along with the tRFs. N: normal BmN cells; V: BmNPV-infected BmN cells. 



stress of viral infection. The levels of the corresponding 
tRNAs were not increased in the infected cells, indicating 
that it is the processing to tRFs that is changed in the 
infected cells. Surprisingly, the ~33nt 5'tRF molecules 
accounted for an overwhelming proportion of the total 
small RNAs bound to BmAgo2 (Figure 2B and Additional 
file 2: Table S2) in the infected cells, especially Bm- 
AspGTC-5'a, which accounted for 49.03% (5,731,905/ 
11,691,441) of the total BmAgo2-associated small RNAs. 



The biogenesis of tRFs under conditions of virus stress pro- 
vides an important piece of evidence for their functional 
role in the cell. A previous report showed that tRFs are 
endogenously associated with Ago2 and are able to guide 
Ago2 to cleave its target RNA [16]. BmAgo2 could thus 
recruit large numbers of tRFs, forming a RNP complex 
similar to RISC and thereby regulating the expression of 
tRF-guided target genes. The precise functions of tRFs 
remain unknown and should be further studied. 
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Transposable-element derived small RNAs and pIRNAs 

To identify small RNAs derived from TEs (TE-derived 
small RNAs) associated with BmAgo2, the unannotated 
reads were mapped to the Bombyx mori TEs. Using this 
method, we obtained 388,788 unique small RNA reads 
derived from the TEs, which accounted for 4778% of 
the total unique reads (388,788/813,702); 56.48% of 
these mapped to multiple locations (219,578/388,788). 
This finding suggests that the population of TE-derived 
small RNAs is significantly higher than that of the small 
RNAs derived from miRNAs, tRNAs, rRNAs, snoRNAs 
and snRNAs in the BmAgo2-associated small RNA 
population. The length distribution ranged from 18nt to 
44nt and peaked at 20nt, with a secondary peak at 21nt, 
implying that the main siRNA proportion comprised 
TE-derived small RNAs (Additional file 2: Table S3). 
TEs are currently regarded as one of the principal forces 
driving genome diversity and evolution but must be 
under appropriate control to maintain genome integrity 
[67,68]. In Bombyx mori, the repetitive sequences are 
estimated to account for 43.6% of the total sequence 
[62,63], implying that active functions and regulation of 
TEs must be in action. We further investigated the TE- 
derived small RNAs restricted to particular classes of 
transposons. Of the obtained 2320 TEs, 2080 could 
generate small RNAs at an abundance between 1 and 
24,484. The majority of these small RNAs were derived 
mainly from LTR and LINE retrotransposons, similar to 
the small RNAs associated with sjAgo2 [18]. The top 25 
TE-derived small RNA clusters and their mapping reads 
are summarized in Table 1 and Additional file 2: Table S4. 
The top 25 TEs consisted of 11 LTR and 9 LINE TEs; 
the only ClassII TE was BmpiggyBac-MER85, a widely 
studied transposon in insects. A previous report showed 
that LINE bml645 strongly and preferentially generated 
antisense small RNAs [18]. We further observed that 
BmAgo2 could associate with small RNAs derived from 
bml645 including 12539 antisense small RNAs and 
1732 sense small RNAs. The other TEs preferentially 
produced antisense small RNAs, including SARTl, 
bm447 and bml770 (antisense/sense > 2). Only Gypsy- 
20_BM-I preferentially produced sense small RNAs 
(sense/antisense >2). 

The piRNAs are generated mainly from germline cells 
and are specifically loaded onto germline-specific 
Argonaute proteins-PIWI proteins [70]. Strangely, using 
piRNApredictor [19] and an in-house developed Perl 
script, a total of 54,033 reads were identified as piRNA 
candidates among the TE-derived small RNA reads 
(Table 1), which accounted for 6.64% of the total unique 
reads. Ping-pong pairs were defined as precise 10-nt 
overlaps between the sense and antisense TE-derived 
piRNAs and were a canonical feature of the piRNAs. 
We also identified 11,337 Ping-pong pairs from the 



54,033 reads, with 5'-U and lOA pairing formed by 
35,630 sense and antisense piRNAs. These results indi- 
cate that a single lOA-piRNA might correspond to more 
than one 5'U-piRNA and vice versa (Additional file 2: 
Table S5). The above results suggest that BmAgo2 may 
function to maintain genome stability by suppressing 
the activities of transposons. 

Known and novel mIRNAs associated with the BmAgo2 
protein 

To function, miRNA must be assembled into RISC [1,71]. 
The Ago protein is the core component of RISC. In recent 
years, the isolation of miRNAs and miRNA targets through 
the immunoprecipitation of Ago protein complexes has 
become a well-characterized method to experimentally 
identify miRNA targets [72,73]. Therefore, large numbers 
of miRNAs should be contained in the BmAgo2-associated 
small RNAs. Of the 567 Bombyx mori miRNAs recorded 
in the miRBase (Release 18, http://www.mirbase.org/), 
270 were identified in the BmAgo2-associated small 
RNA library, including 202 miRNAs, 66 miRNA*s, 1 
miRNA-5p and 1 miRNA-3p species (Additional file 2: 
Table S6). Some of these miRNAs, such as bmo-bantam, 
bmo-miR'184y bmo-miR-2766, bmo-miR-11 and bmo- 
miR-100, exhibited a clear preference for binding to 
BmAgo2 with a high read copy of greater than 5000. 

By mapping all reads to the precursors of known 
miRNAs, we observed that for most known miRNA pre- 
cursors, the matched reads centered around mature miRNA 
sequences and some miRNA* sequences (Additional file 1: 
Figure S4A). However, some reads showed a smear distri- 
bution for certain miRNA precursors (Additional file 1: 
Figure S4B), which was also observed in Arabidopsis Agol 
and Ago4 associated miRNA precursors. These miRNAs 
were not regarded as bona fide miRNAs, but rather 
as siRNAs whose "precursors" also exhibited hairpin 
structures [19]. 

Sequence variants could be involved in the editing 
events, SNPs and their variations. Mature miRNA variants 
are also known as isomiRs [74] and are usually found by 
small RNA deep sequencing [19,75]. A total of 44 different 
isomiRs, each with an abundance of at least 0.02% and 
representing 28 unique miRNA and 1 unique miRNA* 
species, were identified. The sequence length variations 
were not considered in this analysis, and all 44 identified 
isomiR species were represented by different RNA editing 
sites or unannotated SNPs (Additional file 2: Table S7). 
Interestingly, for 8 of these miRNA species, the identified 
responding isomiRs exhibited an abundance of 100%, 
suggesting that the BmAgo2 protein had a clear preference 
for associating with the isomiRs of these miRNAs rather 
than with themselves. 

To identify novel miRNAs, Mireap software was used 
to query the unannotated and other small RNA reads. A 
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Table 1 Top 25 TE-derived small RNAs associated with BmAgo2 



Name of TEs 


TE types 


Length 




piRNAs 


Other small RNAs 








Sense 


Antisense 


Sense 


Antisense 


Gypsy-9_BM(Aquila)-l 


LTR/ Gypsy 


7716 


621(2803) 


692(5297) 


12022(57722)** 


1 1 149(52894) 


bnnl866 


LTR/Pao 


6465 


576(1668) 


613(1988) 


6921(27652) 


7048(28122) 


bnnl645 


LINE/Rl 


13595 


42(54) 


633 (4794) 


1690(3400) 


1 1 906(220080) 


bml456 


LTR/Pao 


8702 


570(2359) 


671(2944) 


4135(22334) 


3842(20498) 


L09635 


LTR/Pao 


4791 


376(1 123) 


408(1644) 


4043(23398) 


4177(18148) 


AB480244.1(BMC1) * 


LINE/Jockey 


25948 


474(1757) 


359(821) 


3598(26788) 


2713(8313) 


AB480245.1(Tanna) * 


non-LTR 


34087 


322(631) 


329(1492) 


1832(6701) 


1 966(4254) 


Yamato 


LTR/Pao 


6400 


207(608) 


299(1628) 


1705(6094) 


2087(9633) 


AB480243.1(Taguchi) * 


LINE/Rl 


32792 


255(738) 


302(841) 


1587(7621) 


1570(5069) 


BnnpiggyBac-MER85 


Classll^:TIR/piggyBac 


3504 


305(1379) 


258(938) 


1 684(9442) 


1456(6983) 


bm228 


LTR/Un known 


8077 


205(473) 


261(786) 


1423(3304) 


1 760(5449) 


BmorL326.1597 


LINE/Unknown 


576 


142(760) 


1 48(800) 


1218(6510) 


2022(17937) 


Gypsy-20_BM-I 


LTR/ Gypsy 


4469 


436(1 125) 


261(3602) 


1974(7179) 


807(41 94) 


AB480240.1(Bml) * 


SINE/SINE 


14359 


337(1093) 


257(1 114) 


1366(3521) 


1466(5126) 


Moriya-I 


LTR/ Bel 


6871 


196(412) 


220 (947) 


1180(3113) 


1350(5170) 


AB480234.1(Bml) * 


SINE/SINE 


18645 


1 69(386) 


1 83(659) 


1095(5477) 


1 246(4397) 


Kabuki 


LTR/Un known 


5342 


221(668) 


240(2406) 


1 1 34(2787) 


1002(2610) 


C A DTI 
oAK 1 1 


LINL/K 1 


D/Uz 


77(1 97) 


1 Zd[Z0U) 


807(1468) 


1575(4629) 


AB480242.1(BMC1) * 


LINE/Jockey 


20035 


223(1024) 


1 76(450) 


1125(4006) 


856(3017) 


bml222 


Unknown/Unknown 


1468 


126(220) 


160 (495) 


1221(3003) 


814(2988) 


bml54 


LINE/Rl 


9083 


104(304) 


73(171) 


1224(4407) 


884(1988) 


AB480235.1 (Kendo) * 


LINE/CRl 


32287 


196(473) 


143(368) 


858(3345) 


622(1317) 


bm447 


LTR/Pao 


7110 


104(140) 


152 (384) 


333(724) 


1197(5351) 


BMC1 


LINE/Jockey 


5091 


85(152) 


120(344) 


765(1623) 


815(2924) 


bml770 


LTR/Pao 


6669 


96(159) 


162(756) 


301(646) 


1089(5433) 



*These contigs from the W chromosome contained TEs and generated a large amount of small RNAs around the TEs, which might be related with the 
transcription of the TEs. 

**Total reads of TE-derived small RNAs are shown in parenthesis. 

^TEs have been divided into two classes based on their transposition intermediate. Class I TEs use an RNA-mediated mode of transposition and are further divided 
into two subclasses, including LTR and non-LTR retrotransposons. Non-LTR retrotransposons are also called LINEs (long interspersed nuclear elements) and SINEs 
(short interspersed nuclear elements). Class II TEs use a DNA-mediated mode of "cut and paste" transposition. In Bombyx mori, many TEs have been identified, 
such as bm447, bm1770, bm1645 and SART1 [69]. 



total of 301 novel miRNA species with a read copy num- 
ber of no less than 5 were identified (Additional file 2: 
Table S6). The precursors of these predicted novel 
miRNAs exhibited the typical stem-loop structure. To 
validate the candidate novel miRNAs, stem-loop qPCR 
assays were developed. The results showed that 7 of the 
8 novel candidate miRNAs selected were detectable in 
the BmN cells, and their expression differed between the 
normal and the BmNPV-infected BmN cells (Additional 
file 1: Figure S5). A total of 301 novel miRNAs were 
identified from among the 892 unique reads. Of the 892 
reads, 356 were derived from unannotated reads; however, 
the remaining 536 were derived from the annotated TE- 
derived small RNAs (469), piRNAs(60) and snoRNA- 
derived small RNAs(7) described above. The TEs and 
snoRNAs could generate functional miRNAs, and this 
phenomenon is regarded as one of the evolutionary 



mechanisms for miRNA formation [76-82]. Previously, 
182 TE-derived miRNAs were also identified [18]. Here, 
our results indicate that Bombyx TEs could generate 
miRNAs associated with the BmAgo2 protein. In addition, 
the association between the RISC effector complex and 
these 301 novel candidate miRNAs also showed the 
biological feature of miRNAs and further validated the 
prediction. 

Viral miRNAs can inhibit both viral and host transcripts, 
and usually help in maintaining the persistent infection of 
a virus or modulating the host antiviral response [83]. 
BmNPV-encoded miRNAs were first identified by J. Singh 
et al [84], and one of these, hmnpv-miR-1, was shown to 
suppress host miRNA biogenesis by regulating the 
exportin-5 cofactor Ran [85]. The 4 miRNAs that were 
first identified did not exist among the 218 BmNPV- 
encoded small RNAs. However, nine novel potential viral 
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miRNAs were identified from among the 218 BmNPV- 
encoded small RNAs (data not shown) using mirExplorer 
software [86], and these should be further validated. 

Other novel small RNAs derived from rRNAs, snoRNAs 
and snRNAs 

Many groups have reported finding novel small RNAs 
derived from housekeeping non-coding RNAs, such as 
rRNAs, snoRNAs and snRNAs [16,87,88]. In addition to 
the tRFs and miRNAs, we found a large number of 
small RNA reads derived from rRNAs, snoRNAs and 
snRNAs, implying that various types of small ncRNAs 
associate with BmAgo2 (see GEO accession numbers 
GSE41841). 

A total of 746 reads with copy numbers less than 50 
mapped to the 5' and 3' termini of 5.8S rRNA. Using 
overlapping probes, the small RNA species derived from 
the 5' termini of 5.8S rRNA were further identified by 
Northern blotting, and their size was found to be approxi- 
mately 50nt. These small RNAs exhibited a slightly lower 
expression level in the virus-infected BmN cells than in 
the normal BmN cells (Figure 4). However, we failed to 



detect small RNAs originating from other locations of 
5.8S rRNA by Northern blotting. 

Li et al previously identified 141 snoRNAs from the 50- 
500nt portion of the total RNA of Bombyx mori using 
RNomics and comparative genomics methods [41]. Of the 
141 snoRNAs, 133 were found to generate small RNAs, 
mainly from 5' and 3' termini. A total of 4205 unique 
reads showed a 5 ' and 3 ' terminal distribution for the 133 
snoRNA sequences (Additional file 2: Table S8). SnoRNAs 
usually function as guide RNAs to direct sequence- 
specific modifications of the targeted RNA. Recent work 
has shown that some miRNA precursors and snoRNAs 
share the same features, and snoRNAs can also give rise 
to small RNAs with miRNA-like functions, suggesting that 
these classes of RNAs may have an evolutionary relation- 
ship [80-82]. Like miRNAs, snoRNA-derived small RNAs 
are also functionally associated with the Argonaute pro- 
tein [88,89]. Using mirExplorer software [86], we found 
that 24 of the 133 snoRNAs exhibited perfect hairpin 
structures similar to those of pre-miRNAs (Additional 
file 2: Table S9); these consisted of 607 unique reads 
with a total read copy number of 107,950 (Additional 




- 5.8S rRNA 
(M36302.1) 



B 



probe p1:GATCCACCGTCCAGGGTAATG 
probe p2:CCGCGAGCCAAGTGATCCACCG 
probe p3:CTTCATCGACCCGCGAGCCAAG 
probe p4:TATGACGCGCAGTTAACTGCGT 




— ^ 5.8S rRNA 
small RNA 



PI P2 P3 P4 

Figure 4 Identification of small RNAs derived from the 5.8S rRNA of Bombyx mori. (A) Mapping of reads to the 5.8S rRNA. A total of 746 
reads with read copy numbers of no less than 50 were mapped mainly to the 5' and 3' termini of the 5.8S rRNA. The designed probes also 
appeared on the corresponding locations. (B) The sequences of the probes used to identify the small RNAs. The pi, p2 and p3 were overlapped 
on the 5.8S rRNA sequence. (C) The small RNAs derived from 5.8S rRNA were identified by Northern blotting. Their size was approximately 50nt. 
The small RNAs originating from the other locations on the 5.8S rRNA were not identified. N: normal BmN cells; V: BmNPV-infected BmN cells. 
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file 2: Table S8). The numerous snoRNA-derived mole- 
cules associated with BmAgo2 will require further study. 

Discussion 

Ago proteins are the key component of the siRNA and 
miRNA pathway and are indispensable binding proteins 
for the function of many other small RNAs. For ex- 
ample, tRF must bind Ago2 to direct the cleavage of 
target RNAs, which contain a perfect complementary 
sequence to the tRF [16]. Using BmAgo2 polyclonal 
antibody, endogenous BmAgo2 was identified in BmN 
cells (Additional file 1: Figure SI). The Ago2 protein 
usually displays a Mg^^-dependent endonuclease activ- 
ity, and Mg^^ provides a favorable contribution to the 
binding of the DNA-RNA to the Ago2 protein [90]. Our 
work demonstrates that recombinant BmAgo2 could 
exhibit an endonuclease activity and cleave the IP RNA 
products when EDTA is removed from the lysis buffer 
(Additional file 1: Figure S6), suggesting recombinant 
HIS -tagged BmAgo2 might exhibit the same physio- 
logical activity as endogenous BmAgo2. 

Using a HIS monoclonal antibody, we successfully iso- 
lated BmAgo2-associated small RNAs by immunopre- 
cipitation from the BmN cells. We used an RIP-seq 
method to obtain a total of 817,302 BmAgo2-associated 
unique small RNAs from the ovary-derived BmN cells, 
91,092 of which we could not match to the reference 
genome of the male silkworm. The unmatched reads 
might be derived from the W chromosome or might 
map to potential gaps in the sequenced female genome. 
To date, many species of small RNAs have been identified 
in Bombyx morU including miRNAs [24-27,29,30,91], 
piRNAs [35,37,38], TE-derived small RNAs [18,34] and 
snoRNA/snRNAs of between 50nt and SOOnt [41]. In the 
present work, we identified BmAgo2-associated small 
RNAs between 18nt and 50nt in length, including the first 
identified tRFs, rRNA-, snoRNA- and snRNA-derived 
small RNAs in this species, as well as previously identified 
miRNAs, piRNAs and TE-derived small RNAs. Import- 
antly, the BmAgo2 protein was found to recruit tRFs at an 
extremely high abundance (Figure 2B), suggesting that 
tRFs could function by binding the BmAgo2 protein. 

We further employed Northern blotting to identify tRFs, 
and the results showed two obvious bands responding to 
the tRNAs and tRFs, respectively, suggesting that tRFs are 
a stable and novel class of small RNAs (Figure 3). Previous 
studies have shown that the ANG nuclease can generate 
30-50nt long sitRNAs under conditions of stress from the 
anti-codon loops, and it has been reported to be a tRNA- 
specific enzyme [92-94]. The Bombyx 5' and 3 'tRFs iden- 
tified here by Northern blotting were approximately ~33nt 
and ~40nt long, respectively. Northern blotting with over- 
lapping probes further confirmed that the ~40nt 3 'tRFs 
were generated by cleaving the anti-codon loops of the 



72nt tRNAs (Additional file 1: Figure S7). These 5' and 3' 
tRFs might be sitRNAs that are produced by the ANG nu- 
clease. Our deep sequencing results further suggested that 
the ~33nt 5 'tRFs and the ~40nt 3 'tRFs could be further 
cleaved on their D and Tv(/C loops by a certain nuclease, 
thereby generating the ~18nt 5' truncated tRFs and ~21nt 
3 'trailer tRFs (Additional file 1: Figure S3), respectively. 
However, we failed to detect these BmAgo2-associated 
tRFs by Northern blotting from among the total RNAs of 
the normal or BmNPV-infected BmN cells, possible due 
to the low abundance of these tRF species. ANG and other 
RNases may be involved in generating these small tRF spe- 
cies, which can cleave tRNAs to produce ~20nt terminal 
tRFs [16]. Furthermore, a previous study also showed that 
the major determinant in the processing of terminal small 
tRFs by ANG and possibly by related RNases is the Tv(/C 
loop within the tRNA. The ANG tRF cleavage site is 
known to be located between the first and second base 
after the GTTC motif on the Tv|/C loop [16]. Our deep se- 
quencing data show that a conserved GTTC domain exists 
at the termini of most 21nt 3 'trailer tRFs, suggesting there 
is a potential tRF cleavage site between the GTTC cytosine 
and the first base after the GTTC motif on Tv|/C loop 
(Additional file 1: Figure S3). The 3 'trailer tRFs could 
not be detected by Northern blotting as they had a low 
total RNA abundance; however, the 3' trailer tRFs were 
significantly more abundant than the 3 'tRFs in the 
BmAgo2-associated total RNA. In the sequenced 
BmAgo2-associated small RNAs, the read copies of the 
3 'trailer tRFs were far greater than those of the 3 'tRFs 
(Additional file 2: Table S2), suggesting that the 3 'trailer 
tRFs preferentially associated with BmAgo2 for their 
functions, which might be similar to those of the siRNA 
pathway. For the tRNAs with a variable arm, the D loop 
could be cleaved at the GGCCGAGCGG cleavage site, 
generating the ~18nt 5 'truncated tRFs (Additional file 1: 
Figure S3D). This result implies a possible relationship 
between the D loop cleavage and the variable arm during 
the biogenesis of 5 'truncated tRFs. 

The tRF biogenesis that occurs under conditions of viral 
stress provides an important piece of evidence about the 
functional role of tRFs. The present work shows that when 
the BmN cells were infected with BmNPV, many tRFs 
were generated and recruited onto the BmAgo2 protein. 
Therefore, these tRFs were enriched among the BmAgo2- 
associated small RNAs. The majority of these 5 'tRFs iden- 
tified by Northern blotting were expressed only in 
BmNPV-infected cells (Figure 3A), whereas the other tRF 
species identified by Northern blotting were usually 
expressed in both the normal and the virus-infected cells 
(Figure 3B), suggesting the main tRF species generated by 
viral stress may be the 5 'tRFs. The asymmetric expression 
of 5 'tRFs implies a relationship between the 5 'tRFs and 
BmNPV infection in the BmN cells. Furthermore, 5 'tRFs 
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accounted for an overwhelming proportion of the 
BmAgo2-associated small RNAs. This was especially true 
for Bm-AspGTC-5'a, which accounted for 59.58% of the 
total BmAgo2-associated small RNAs. These observations 
suggest that the 5'tRFs are likely to bind the BmAgo2 pro- 
tein and thus influence the BmNPV infection. 

The piRNA system is thought to provide a germline 
defense against TE activity. A simflar function was also 
recently noted in the siRNA pathway [51,68,95-97]. 
Drosophfla utilizes two small-RNA systems to restrict 
transposon activity in its germline (mostly via piRNAs) 
and soma (mostly via siRNAs); in this system, Ago2 
deficiency results in increased levels of transposon 
transcripts [95]. In mouse oocytes, the loss of Ago2 can 
result in decreased levels of siRNAs and increased levels 
of retrotransposons [97]. In S. japonicum, sjAgo2 
strongly associated with TE-derived siRNAs and is be- 
lieved to regulate transposons [18]. These observations 
indicate that the Ago2 protein might be involved in 
regulating transposons at the transcriptional level in 
germline cells. The BmN cell line was obtained from 
the Bombyx ovary. TE-derived small RNAs accounted 
for almost half of the total unique BmAgo2-associated 
small RNAs. The peak of their length distribution was 
the typical siRNA length. Therefore, we suggested that 
BmAgo2 could execute important functions in the main- 
tenance of genome stability in germline cells through an 
RNAi-like pathway. 

piRNApredictor is a software program for piRNA pre- 
diction that uses a novel 'dynamic' algorithm with a pre- 
cision of over 90% and a sensitivity of over 60% [19]. 
Using piRNApredictor supplemented with a perl script 
for screening Ping-pong pairs, 54,033 >25nt piRNA can- 
didates and 11,337 typical Ping-pong pairs were identified 
from among the TE-derived small RNAs associated with 
BmAgo2. In addition, we further identified 56,131 >25nt 
piRNA candidates and 11,898 Ping-pong pairs with 17,533 
sense piRNAs and 18,239 antisense piRNAs from 
the remaining unannotated reads (Additional file 2: 
Table SIO). These piRNA candidates could not be 
mapped to the identified TEs and could be derived from 
some heterochromatic regions, other repeat sequences 
and unknown transposons in the Bombyx genome. piRNA 
is thought to defend the host genome against transposons 
in Bombyx mori [34]. Recently, an increasing number of 
piRNAs have been identified in Bombyx mori [18,35,37], 
and their biogenesis has been further studied, providing 
evidence for their "Ping-pong Circle" and biosynthesis 
mechanisms [37-40]. We further compared our data with 
the previously identified piRNAs [18,35,37]. The results 
showed that our data has 0.83% (909), 14.95% (16,474), 
30.62% (33,733), 36.54% (40,254), 17.55% (19,335) and 
17.72% (19,517) overlap (mismatch not allowed) with the 
TE-, BmN4-, Siwi-, BmAgo3-, ovary- and testis-derived 



piRNAs, respectively. This overlap serves as an indicator 
of the high reliability of our dataset. We first identified 
piRNAs associated with the BmAgo2 protein, which also 
contains PAZ and PIWI domains. BmAgo2 might bind to 
piRNA-like RNA molecules through the PAZ domain, the 
biogenesis of which could be similar with that of piRNAs. 
The observations above suggest that BmAgo2 can associ- 
ate with the 20 22nt TE-derived siRNAs and the >25nt 
piRNA-like small RNAs during the regulation of trans- 
poson expression. 

We must emphasize that many more reads did not 
align to previously known ncRNAs, such as miRNAs in 
miRBase, tRNAs, rRNAs, piRNAs, snoRNAs, snRNAs 
and other reported Bombyx mori small RNAs [41], in- 
dicating that many small RNAs remain to be discov- 
ered in Bombyx mori. The unannotated small RNAs 
associated with BmAgo2 usually exhibit a higher abun- 
dance than the miRNAs. Only 49 miRNAs had a read 
copy of at least 100 (Additional file 2: Table S6), but 
848 unannotated small RNAs had a read copy of at 
least 100. 

Conclusions 

We successfully overexpressed BmAgo2 under the control 
of the iel promoter in BmN cells using a recombinant 
BmNPV virus and used an RIP-seq approach to identify 
the associated small non-coding RNAs (18-50nt). Various 
types of small ncRNAs were found to be associated with 
BmAgo2, including tRNA-, TE-, rRNA-, snoRNA- and 
snRNA-derived small RNAs, as well as miRNAs and 
piRNAs. The tRFs were the most abundant class of 
BmAgo2-associated small RNAs. However, Northern blot- 
ting analysis showed that many of the tRFs were only 
expressed or significantly up-regulated in the BmNPV 
infected cells, implying that they act by binding BmAgo2 
during the response of Bombyx mori to BmNPV. We also 
demonstrated that a significant fraction of the small RNAs 
associated with BmAgo2 are derived from TEs, many of 
which were identified as piRNA candidates. This finding 
suggests that they may function to maintain genome 
stability through the suppression of transposon activity. 
Finally, we also identified a large number of miRNAs 
associated with BmAgo2, including 301 novel miRNAs. 
Taken together, the above findings provide new clues 
for the future functional study of small RNAs during 
insect development and evolution. 

Methods 

Cell line 

The BmN ceU line, maintained in our laboratory and de- 
rived from the silkworm ovary, was cultured in Sf-900 II 
SFM medium (Gibco BRL) supplemented with 10% (v/v) 
fetal bovine serum (Gibco BRL) at 27°C. 
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Expression of BmAgo2 using the Bacmid expression 
system and Western blotting analysis 

BmAgo2 was expressed using the Bacmid system under 
the control of the iel promoter enhanced with hrS of 
Autographa Calif ornica nucleopolyhedro virus (AcNPV). 
The recombinant Bacmid vector harboring the BmAgo2 
gene was previously constructed by us and was preserved 
in our laboratory [61]. In this Bacmid, the immediate early 
iel promoter and the hrS enhancer were introduced and 
used to drive the expression of BmAgo2. BmN cells were 
infected with the constructed recombinant virus at an 
MOI of 10 and cultured at 27°C. After 17 h, the cells were 
harvested and washed twice with PBS (pH 7.4). The cells 
were then lysed in cold lysis buffer (50 mM Tris, pH 8.0, 
150 mM NaCl, 5 mM EDTA, 0.5% Nonidet P-40) with 
freshly added 1 mM DTT for 30 min. The expressed 
HIS-BmAgo2 protein was analyzed using the previously 
described method of Western blotting with the mouse 
Anti-Hise monoclonal antibody (Roche) (1:2000) [98]. 
Briefly, cell extracts with over-expressed BmAgo2 were 
separated on an SDS-PAGE gel and transferred to a PVDF 
membrane (Millipore). The membrane was blocked with 
5% BSA for 90 min at room temperature. Anti-Hise anti- 
body was used for the detection of HIS-BmAgo2. The 
HRP-conjugated goat mouse IgG (1:10,000) was used as 
the secondary antibody, and the signal was detected using 
the DAB kit (Gersion). 



Co-immunoprecipitation 

Cells infected with recombinant virus were lysed in an 
ice-cold lysis buffer (50 mM Tris, pH 8.0, 150 mM NaCl, 
5 mM EDTA, 0.5% Nonidet P-40) with freshly added 
1 mM DTT, 200 U/ml RNase inhibitor (Promega) and a 
protease inhibitor cocktail (Roche) for 30 min at -80°C. 
After dissolving at room temperature, centrifugation was 
performed at 12000 rpm at 4°C for 10 min. The super- 
natant was first mixed with 40 [A of protein G magnetic 
beads (Millipore) and incubated at 4°C for 3 h with gen- 
tle agitation (Mock). After immobilizing the magnetic 
beads with a magnet, the supernatant was removed and 
subsequently mixed with 5 [A of Anti-Hisg antibody or 
normal mouse IgG for 3 h at 4°C with gentle agitation. 
Then, 40 [A of the protein G magnetic beads were added 
and continually incubated for 3 h at 4°C. The immunopre- 
cipitate was washed five times in lysis buffer and separated 
into two portions, one for the RNA isolation as a sample 
for the following experiment and another for the Western 
blotting to identify the immunoprecipitation of BmAgo2. 
Normal Mouse IgG was used as a negative control in the 
immunoprecipitation procedure. The eluted proteins were 
separated on SDS-PAGE gel and identified by Western 
blotting with an Anti-Hise antibody (1:2000) and the ECL 
Plus Western Blotting Detection System (GE Healthcare). 



RNA extraction 

The total RNAs were extracted from the immunopre- 
cipitate of BmAgo2 using TRIzol reagent (Invitrogen) 
according to the manufacturer s protocol, the concen- 
tration and absorbance of which were determined using 
a Nanodrop ND-1000 spectrophotometer. The total 
RNAs were subsequently loaded onto a 15% denaturing 
polyacrylamide gel (PAGE) containing 8 M urea. They 
were electrophoresed and then stained with 1 [ig/[A 
ethidium bromide. 

RNA-seq 

A small RNA library was constructed using the total 
RNAs extracted above. Briefly, the RNAs were size- 
fractionated on a 15% polyacrylamide gel, and the 18- 
50nt fraction was collected. The collected small RNAs 
were ligated with 5 ' and 3 ' Illumina adaptors and sub- 
sequently used as a template to synthesize first-strand 
cDNA. The cDNA was amplified by PCR with the 
Illumina small RNA primer set and sequenced on the 
Illumina HiSeq 2000 at BGI (Beijing Genomics Institute, 
Shenzhen, China). We failed to construct a small RNA 
library for the negative control due to the low concentra- 
tion of total RNAs extracted from the negative control of 
the immunoprecipitate. 

Small RNA sequence analysis 

Adaptors, low-quality tags and contaminants were re- 
moved from raw reads produced by Solexa sequencing 
to produce the clean total reads. All identical reads were 
counted and merged as unique reads. The length distri- 
bution of the clean total reads and unique reads were 
then summarized. These unique reads with read counts 
were mapped to the Bombyx mori [62,63] and BmNPV 
genome [99] using Bowtie or SOAP2 software with 0 or 
2nt mismatch [100,101]. Afterwards, the clean reads 
that perfectly mapped the genome were annotated into 
different categories according to the standard bioinfor- 
matics pipeline designed by GBI. Briefly, by comparing 
the clean reads with those in existing databases and 
picking out the overlap on genome location between our 
data and the databases using BLAST or other programs 
developed by GBI, the small RNA reads were annotated 
into different categories. These categories included 
miRNAs as well as TE-, tRNA-, rRNA-, snoRNA- and 
snRNA-derived small RNAs. The remaining reads were 
then used to predict the novel miRNAs and to base edit 
the potential known miRNAs. The >25nt TE-derived small 
RNAs and remaining unannotated reads were selected for 
probing piRNAs using the piRNApredictor program [19] 
and Perl script developed by us. The reads that perfectly 
mapped to the BmNPV genome were used to predict the 
BmNPV-encoded miRNAs using mirExplorer software 
[86]. Finally, the remaining and unmapped reads were 
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labeled as unannotated small RNAs. The whole process is 
shown in Additional file 1: Figure S8. 

Datasets downloaded for the bioinformatics analysis 
above are listed as follows: Bombyx miRNAs was derived 
from miRbase 18.0 (http://www.mirbase.org). Silkworm 
ncRNAs used for the annotations, including the tRNA, 
rRNA, scRNA, snoRNA and snRNA annotations, were 
obtained from NCBI and Rfam (http://www.ncbi.nlm. 
nih.gov and http://rfam.sanger.ac.uk); 1476 annotated 
silkworm TE sequences were retrieved from NCBI, 
MAS DNA Bank (http://www.dna.affrc.go.jp) and the 
BmTEdb database (http://gene.cqu.edu.cn/BmTEdb) and 
1668 Re AS clones were obtained from the BmTE Library 
(http://sgp.dna.affrc.go.jp/pubdata/genomicsequences.html). 
The obtained silkworm TE sequences were merged 
using CD-HIT software [102]. Some contigs containing 
TEs from the W chromosome were also included in our 
TE datasets to probe TE-derived small RNAs, which 
might be related to the transcription of the TEs. The 
BmN4-, Siwi-, BmAgoS-, ovary-, testis- and TE-derived 
piRNAs used for the comparative analysis were retrieved 
from DDBJ (http://www.ddbj.nig.ac.jp/) and the Supple- 
mentary materials according to the references [18,35,37]. 

Northern blotting 

The oligonucleotide probes were labeled with the DIG 
Oligonucleotide Tailing Kit (Roche) according to the 
manufacturers instructions. Illustra ProbeQuant G-25 
Micro Columns (GE Healthcare) were used for the puri- 
fication of probes after labeling. 

Northern blotting analysis was performed on the total 
RNA extracted from the normal BmN cells and the 
BmN cells infected with the BmNPV virus. The total 
RNAs were mixed with Gel Loading Buffer II (Ambion), 
heated to 95°C for 3 min and separated on 15% denatur- 
ing TBE-Urea gels. Next, the RNAs were transferred 
onto the charged nylon membranes (Millipore), which 
were then subjected to UV cross-linking, pre-hybridized 
at 42°C for 2 h using the PerfectHyb Plus Hybridization 
Buffer (Sigma) and hybridized with probes overnight at 
37°C. After hybridization, the membranes were rinsed and 
washed four times (twice at 25°C for 15 min using 2x SSC 
and 0.1% SDS, twice at 37°C with 0.5x SSC and 0.1% 
SDS). Signals were detected with a DIG Luminescent 
Detection Kit (Roche) according to the manufacturers 
protocol and were finally exposed to X-OMAT BT film 
(Kodak). 

Stem-loop qRT-PCR 

The novel miRNAs predicted were further identified by 
miRNA-specific stem-loop quantitative real-time PCR 
using the mirVana™ qRT-PCR miRNA Detection kit 
(Ambion) as described previously [103]. The top 8 novel 



miRNAs with the highest read copy number were se- 
lected for qRT-PCR identification. 

Additional files 



Additional file 1: Supplementary figures contain 8 figures. 

Additional file 2: Annotation and analysis of small RNAs associated 
with BmAgo2 contain 10 tables. 



Abbreviations 

ncRNAs: non-coding RNAs; BmAgo2: Bombyx mori argonaute2; tRFs: tRNA- 
derived fragments; miRNAs: microRNAs; siRNAs: small-interfering RNAs; 
piRNAs: piwi-interacting RNAs; IP: Immunoprecipitation; RIP: RNA 
immunoprecipitation; TE: Transposable element; RISC: RNA-induced silencing 
complex; AcNPV: Autographa Californica nucleopolyhedrovirus; 
BmNPV: Bombyx mori nucleopolyhedrovirus; PAGE: Polyacrylamide gel. 

Competing interests 

The authors declare that they have no competing interests. 
Authors' contributions 

ZMN performed the experimental design and bioinformatics analysis and 
drafted the manuscript. FZ performed the IP and the construction of the 
small RNA libraries. DL, ZBL and JC performed the Northern blotting analysis. 
YL and JHS performed the stem-loop qRT-PCR analysis. QS performed the IP. 
WY was responsible for the cell culture. WPZ and CYJ performed the 
expression study of BmAgo2. YHY performed the bioinformatics analysis. JMY 
and YFJ designed the small RNA sequence analysis. YZZ was responsible for 
the experimental design and revised the manuscript. All authors approved 
the final version of the paper. 

Acknowledgements 

This work was supported by financial grants from the National Basic 
Research Program of China (201 2CB1 14600), the Zhejiang Natural Science 
Foundation (LY12C06003), the National High Technology Research and 
Development Program (201 lAA 100603), the National Natural Science 
Foundation of China (31125011) and the "521" talents development program 
of Zhejiang Sci-Tech University. We thank Prof Qingyou Xia for his kind 
discussion about the small RNAs of silkworm and the members of the 
Beijing Computing Center for their assistance with the high-throughput data 
analysis. 

Author details 

^College of Life Sciences, Zhejiang Sci-Tech University, Hanghzou 310018, 
China. ^Zhejiang Economic & Trade Polytechnic, Hangzhou 310018, China. 
^College of Materials and Textile, Zhejiang Sci-Tech University, Hangzhou 
310018, China. ^College of Life Sciences, Zhejiang University, Hangzhou 
310058, China. 

Received: 27 April 2013 Accepted: 26 September 2013 
Published: 28 September 2013 

References 

1. Krol J, Loedige I, Filipowicz W: The widespread regulation of microRNA 
biogenesis, function and decay. Nat Rev Genet 2010, 1 1(9):597-610. 

2. Cai Y, Yu X, Hu S, Yu J: A brief review on the mechanisms of miRNA 
regulation. Genomics Proteomics Bioinform 2009, 7(4):147-154. 

3. Doench JG, Petersen CP, Sharp PA: siRNAs can function as miRNAs. Genes 
Dev 2003, 1 7(4):438-442. 

4. Sharma N, Sahu PP, Puranik S, Prasad M: Recent Advances in plant-virus 
interaction with emphasis on small interfering RNAs (siRNAs). Mol 
Biotechnol 2013, 55(l):63-77. 

5. Gonzalez S, Pisano DG, Serrano M: Mechanistic principles of chromatin 
remodeling guided by siRNAs and miRNAs. Cell Cycle 2008, 

7(1 6):2601 -2608. 

6. Klenov MS, Lavrov SA, Stolyarenko AD, Ryazansky SS, Aravin AA, TuschI T, 
Gvozdev VA: Repeat-associated siRNAs cause chromatin silencing of 
retrotransposons in the Drosophila melanogaster germline. Nucleic Acids 
Res 2007, 35(16):5430-5438. 



Nie et al. BMC Genomics 2013, 14:661 
httpy/www.biomedcentral.com/l 471 -21 64/1 4/661 



Page 14 of 15 



7. Zhang X, Xia J, Lii YE, Barrera-Figueroa BE, Zhou X, Gao S, Lu L, Niu D, Chen 
Z, Leung C, et al: Genome-wide analysis of plant nat-siRNAs reveals 
insights into their distribution, biogenesis and function. Genome Biol 

2012, 13(3):R20. 

8. Smalheiser NR: The search for endogenous siRNAs in the mammalian 

brain. Exp Neurol 2012, 235(2):455-463. 

9. Snead NM, Rossi JJ: Biogenesis and function of endogenous and 
exogenous siRNAs. Wiley Interdiscip Rev RNA 2010, 1(1):1 17-131. 

10. Tushif JS, Zamore PD, Zhang Z: Snapshot: mouse piRNAs, PIWI proteins, 
and the ping-pong cycle. Cell 2009, 1 39(4):830-830. e831. 

11. Tushir JS, Zamore PD, Zhang Z: Snapshot: Fly piRNAs, PIWI proteins, and 
the ping-pong cycle. Cell 2009, 1 39(3):634-634. e631. 

12. Shpiz S, Olovnikov I, Sergeeva A, Lavrov S, Abramov Y, Savitsky M, 
Kalmykova A: Mechanism of the piRNA-mediated silencing of Drosophila 
telomeric retrotransposons. Nucleic Acids Res 201 1, 39(20):8703-871 1. 

13. Castaneda J, Genzor P, Bortvin A: piRNAs, transposon silencing, and 
germline genome integrity. Mutat Res 201 1, 714(1-2):95-104. 

14. O'Donnell KA, Boeke JD: Mighty Piwis defend the germline against 
genome intruders. Cell 2007, 129(1):37-44. 

15. Pillai RS, Chuma S: piRNAs and their involvement in male germline 
development in mice. Dev Growth Differ 2012, 54(1):78-92. 

16. Li Z, Ender C, Meister G, Moore PS, Chang Y, John B: Extensive terminal 
and asymmetric processing of small RNAs from rRNAs, snoRNAs, snRNAs, 
and tRNAs. Nucleic Acids Res 2012, 40(14):6787-6799. 

17. Heyer R, Dorr M, Jellen-Ritter A, Spath B, Babski J, Jaschinski K, Soppa J, 
Marchfelder A: High throughput sequencing reveals a plethora of small 
RNAs including tRNA derived fragments in Haloferax volcanii. RNA Biol 
2012, 9(7):101 1-1018. 

18. Cai P, Piao X, Hou N, Liu S, Wang H, Chen Q: Identification and 
characterization of Argonaute protein, Ago2 and its associated small 
RNAs in Schistosoma japonicum. PLoS Negl Trop Dis 2012, 6(7):e1745. 

19. Zhang Y, Wang X, Kang L: A k-mer scheme to predict piRNAs and 
characterize locust piRNAs. Bioinform 2011, 27(6)771-776. 

20. Sobala A, Hutvagner G: Transfer RNA-derived fragments: origins, 
processing, and functions. Wiley Interdiscip Rev RNA 201 1, 2(6):853-862. 

21. Burroughs AM, Ando Y, de Hoon MJ, Tomaru Y, Suzuki H, Hayashizaki Y, 
Daub CO: Deep-sequencing of human Argonaute-associated small RNAs 
provides insight into miRNA sorting and reveals Argonaute association 
with RNA fragments of diverse origin. RNA Biol 201 1, 8(1 ):1 58-1 77. 

22. Mita K: Genome of a lepidopteran model insect, the silkworm Bombyx 
mori. Seikagaku 2009, 81(5):353-360. 

23. Hamamoto H, Tonoike A, Narushima K, Horie R, Sekimizu K: Silkworm as a 
model animal to evaluate drug candidate toxicity and metabolism. 
Comp Biochem Physiol C Toxicol Pharmacol 2009, 149(3):334-339. 

24. He PA, Nie Z, Chen J, Lv Z, Sheng Q, Zhou S, Gao X, Kong L, Wu X, Jin Y, 
et al: Identification and characteristics of microRNAs from Bombyx mori. 
BMC Genomics 2008, 9:248. 

25. Cai Y, Yu X, Zhou Q, Yu C, Hu H, Liu J, Lin H, Yang J, Zhang B, Cui P, et al: 
Novel microRNAs in silkworm {Bombyx mori). Funct Integr Genomics 2010, 
10(3):405-415. 

26. Yu X, Zhou Q, Li SC, Luo Q, Cai Y, Lin WC, Chen H, Yang Y, Hu S, Yu J: The 
silkworm {Bombyx mori) microRNAs and their expressions in multiple 
developmental stages. PLoS One 2008, 3(8):e2997. 

27. Yu X, Zhou Q, Cai Y, Luo Q, Lin H, Hu S, Yu J: A discovery of novel microRNAs 
in the silkworm {Bombyx mori) genome. Genomics 2009, 94(6):438-444. 

28. Huang Y, Zou Q, Tang SM, Wang LG, Shen XJ: Computational 
identification and characteristics of novel microRNAs from the silkworm 
{Bombyx mori L). Mol Biol Rep 2010, 37(7):3171-3176. 

29. Liu S, Li D, Li Q, Zhao P, Xiang Z, Xia Q: MicroRNAs of Bombyx mori 
identified by Solexa sequencing. BMC Genomics 2010, 1 1:148. 

30. Liu S, Gao S, Zhang D, Yin J, Xiang Z, Xia Q: MicroRNAs show diverse and 
dynamic expression patterns in multiple tissues of Bombyx mori. BMC 
Genomics 2010, 11:85. 

31 . Zhang Y, Zhou X, Ge X, Jiang J, Li M, Jia S, Yang X, Kan Y, Miao X, Zhao G, 
et al: Insect-Specific microRNA involved in the development of the 
silkworm bombyx Mori. PLoS One 2009, 4(3):e4677. 

32. Liu S, Xia Q, Zhao P, Cheng T, Hong K, Xiang Z: Characterization and 
expression patterns of let-7 microRNA in the silkworm {Bombyx mori). 
BMC Dev Biol 2007, 7:88. 

33. Jagadeeswaran G, Zheng Y, Sumathipala N, Jiang H, Arrese EL, Soulages JL, 
Zhang W, Sunkar R: Deep sequencing of small RNA libraries reveals 



dynamic regulation of conserved and novel microRNAs and microRNA- 
stars during silkworm development. BMC Genomics 2010, 1 1:52. 

34. Kawaoka S, Hayashi N, Katsuma S, Kishino H, Kohara Y, Mita K, Shimada T: 
Bombyx small RNAs: genomic defense system against transposons in the 
silkworm, Bombyx mori. Insect Biochem Mol Biol 2008, 38(1 2):1 058-1 065. 

35. Kawaoka S, Hayashi N, Suzuki Y, Abe H, Sugano S, Tormari Y, Shimada T, 
Katsuma S: The Bombyx ovary-derived cell line endogenously expresses 
PlWI/PIWI-interacting RNA complexes. Rna 2009, 15(7):1 258-1 264. 

36. Kawaoka S, Kadota K, Arai Y, Suzuki Y, Fujii T, Abe H, Yasukochi Y, Mita K, 
Sugano S, Shimizu K, et al: The silkworm W chromosome is a source of 
female-enriched piRNAs. Rna 2011, 1 7(1 2):21 44-21 51. 

37. Kawaoka S, Arai Y, Kadota K, Suzuki Y, Hara K, Sugano S, Shimizu K, Tomari 
Y, Shimada T, Katsuma S: Zygotic amplification of secondary piRNAs 
during silkworm embryogenesis. Rna 201 1, 17(7):1 401 -1407. 

38. Kawaoka S, Mitsutake H, Kiuchi T, Kobayashi M, Yoshikawa M, Suzuki Y, Sugano 
S, Shimada T, Kobayashi J, Tomari Y, et al: A role for transcription from a 
piRNA cluster in de novo piRNA production. Rna 2012, 18(2):265-273. 

39. Kawaoka S, Hara K, Shoji K, Kobayashi M, Shimada T, Sugano S, Tomari Y, 
Suzuki Y, Katsuma S: The comprehensive epigenome map of piRNA 
clusters. Nucleic Acids Res 2013, 41(3):1581-1590. 

40. Hara K, Fujii T, Suzuki Y, Sugano S, Shimada T, Katsuma S, Kawaoka S: 
Altered expression of testis-specific genes, piRNAs, and transposons in 
the silkworm ovary masculinized by a W chromosome mutation. BMC 
Genomics 2012, 13:119. 

41 . Li D, Wang Y, Zhang K, Jiao Z, Zhu X, Skogerboe G, Guo X, Chinnusamy V, Bi L, 
Huang Y, etal: Experimental RNomics and genomic comparative analysis 
reveal a large group of species-specific small non-message RNAs in the 
silkworm Bombyx mori. Nucleic Acids Res 201 1 , 39(9):3792-3805. 

42. Bartel DP: MicroRNAs: genomics, biogenesis, mechanism, and function. 
Ce//2004, 116(2):28 1-297. 

43. Mallory A, Vaucheret H: Form, function, and regulation of ARGONAUTE 
proteins. Plant Cell 2010, 22(12):3879-3889. 

44. O'Carroll D, Mecklenbrauker I, Das PP, Santana A, Koenig U, Enright AJ, 
Miska EA, Tarakhovsky A: A Slicer-independent role for Argonaute 2 in 
hematopoiesis and the microRNA pathway. Genes Dev 2007, 

21 (16):1 999-2004. 

45. Hammond SM, Boettcher S, Caudy AA, Kobayashi R, Hannon GJ: 
Argonaute2, a link between genetic and biochemical analyses of RNAi. 
Sci Signal 2001, 293(5532):1 146. 

46. Peters L, Meister G: Argonaute proteins: mediators of RNA silencing. 
Mol Cell 2007, 26(5):61 1-623. 

47. Filipowicz W, Bhattacharyya SN, Sonenberg N: Mechanisms of post- 
transcriptional regulation by microRNAs: are the answers in sight? 
Nat Rev Genet 2008, 9(2):102-1 14. 

48. Dueck A, Ziegler C, Eichner A, Berezikov E, Meister G: microRNAs associated 
with the different human Argonaute proteins. Nucleic Acids Res 2012, 
40(19):9850-9862. 

49. Gu W, Claycomb JM, Batista PJ, Mello CC, Conte D: Cloning Argonaute- 
associated small RNAs from Caenorhabditis elegans. Methods Mol Biol 
2011,725:251-280. 

50. Qi Y, Mi S: Purification of Arabidopsis argonaute complexes and 
associated small RNAs. Methods Mol Biol 2010, 592:243-254. 

51 . Kawamura Y, Saito K, Kin T, Ono Y, Asai K, Sunohara T, Okada TN, Siomi MC, 
Siomi H: Drosophila endogenous small RNAs bind to Argonaute 2 in 
somatic cells. Nature 2008, 453(71 96):793-797. 

52. Flores-Jasso CF, Salomon WE, Zamore PD: Rapid and specific purification 
of Argonaute-small RNA complexes from crude cell lysates. Rna 2013, 
19(2):271-279. 

53. Beitzinger M, Peters L, Zhu JY, Kremmer E, Meister G: Identification of 
human microRNA targets from isolated argonaute protein complexes. 
RNA Biol 2007, 4(2):76-84. 

54. Beitzinger M, Meister G: Experimental identification of microRNA targets 
by immunoprecipitation of Argonaute protein complexes. Methods Mol 
Biol 2011, 732:153-167. 

55. Tsukioka H, Takahashi M, Mon H, Okano K, Mita K, Shimada T, Lee JM, 
Kawaguchi Y, Koga K, Kusakabe T: Role of the silkworm argonaute2 
homolog gene in double-strand break repair of extra chromosomal 
DNA. Nucleic acids Res 2006, 34(4):1 092-1 101. 

56. Wang G, Jiang L, Zhu L, Cheng T, Niu W, Yan Y, Xia Q: Characterization of 
Argonaute family members in the silkworm, Bombyx mori. Insect Sci 2013, 
20(1):78-91. 



Nie et al. BMC Genomics 2013, 14:661 
httpy/www.biomedcentral.com/l 471 -21 64/1 4/661 



Page 15 of 15 



57. Song JJ, Smith SK, Hannon GJ, Joshua-Tor L: Crystal structure of Argonaute 
and its implications for RISC slicer activity. Science 2004, 

305(5689):! 434-1 437. 

58. Yuan YR, Pei Y, Ma JB, Kuryavyi V, Zhadina M, Meister G, Chen HY, Dauter Z, 
TuschI T, Patel DJ: Crystal structure of A. aeolicus argonaute, a site- 
specific DNA-guided endoribonuclease, provides insights into RISC- 
mediated mRNA cleavage. Mol Cell 2005, 1 9(3):405-41 9. 

59. Jayaseelan S, Doyle F, Currenti S, Tenenbaum SA: RIP: an mRNA localization 
technique. Methods Mol Biol 201 1, 714:407-422. 

60. Jain R, Devine T, George AD, Chittur SV, Baroni TE, Penalva LO, Tenenbaum 
SA: RIP-Chip analysis: RNA-Binding Protein Immunoprecipitation- 
Microarray (Chip) Profiling. Methods Mol Biol 201 1, 703:247-263. 

61 . Zhou F, Gao Z, Lv Z, Chen J, Hong Y, Yu W, Wang D, Jiang C, Wu X, Zhang 
Y, et al: Construction of the iel-bacmid expression system and its use to 
express EGFP and BmAG02 in BmN cells. AppI Biochem Biotechnol 2013, 
169(8):2237-2247. 

62. Xia Q, Zhou Z, Lu C, Cheng D, Dai F, Li B, Zhao P, Zha X, Cheng T, Chai C, 
et al: A draft sequence for the genome of the domesticated silkworm 
{Bombyx mori). Science 2004, 306(5703):1 937-1 940. 

63. Mita K, Kasahara M, Sasaki S, Nagayasu Y, Yamada T, Kanamori H, Namiki N, 
Kitagawa M, Yamashita H, Yasukochi Y, et al: The genome sequence of 
silkworm, Bombyx mori. DNA Res 2004, 1 1(1):27-35. 

64. Mi S, Cai T, Hu Y, Chen Y, Hodges E, Ni F, Wu L, Li S, Zhou H, Long C, et al: 
Sorting of small RNAs into Arabidopsis argonaute complexes is directed 
by the 5' terminal nucleotide. Cell 2008, 1 33(1 ):1 16-127. 

65. Takeda A, Iwasaki S, Watanabe T, Utsumi M, Watanabe Y: The mechanism 
selecting the guide strand from small RNA duplexes is different among 
argonaute proteins. Plant Cell Physiol 2008, 49(4):493-500. 

66. Pederson T: Regulatory RNAs derived from transfer RNA? Rna 2010, 
1 6(1 0):1 865-1 869. 

67. Biemont C, Vieira C: Genetics: junk DNA as an evolutionary force. Nature 
2006, 443(71 11):521 -524. 

68. Obbard DJ, Finnegan DJ: RNA interference: endogenous siRNAs derived 
from transposable elements. Curr Biol 2008, 18(13):R561-R563. 

69. Osanai-Futahashi M, Suetsugu Y, Mita K, Fujiwara H: Genome-wide 
screening and characterization of transposable elements and their 
distribution analysis in the silkworm, Bombyx mori. Insect Biochem Mol Biol 
2008, 38(1 2):1 046-1 057. 

70. Siomi MC, Miyoshi T, Siomi H: piRNA-mediated silencing in Drosophila 
germlines. Semin Cell Dev Biol 2010, 21 (7):754-759. 

71. Winter J, Jung S, Keller S, Gregory Rl, Diederichs S: Many roads to maturity: 
microRNA biogenesis pathways and their regulation. Nat Cell Biol 2009, 

1 1 (3):228-234. 

72. Hong X, Hammell M, Ambros V, Cohen SM: Immunopurification of Agol 
miRNPs selects for a distinct class of microRNA targets. Proc Natl Acad Sci 
USA 2009, 106(35):1 5085-1 5090. 

73. Matkovich SJ, Van Booven DJ, Eschenbacher WH, Dorn GW 2nd: RISC RNA 
sequencing for context-specific identification of in vivo microRNA 
targets. Ore Res 201 1, 1 08(1 ):1 8-26. 

74. Landgraf P, Rusu M, Sheridan R, Sewer A, lovino N, Aravin A, Pfeffer S, Rice 
A, Kamphorst AO, Landthaler M, et al: A mammalian microRNA expression 
atlas based on small RNA library sequencing. Cell 2007, 129(7):1401-1414. 

75. Voellenkle C, Rooij J, Guffanti A, Brini E, Fasanaro P, Isaia E, Croft L, David M, 
Capogrossi MC, Moles A, et al: Deep-sequencing of endothelial cells 
exposed to hypoxia reveals the complexity of known and novel 
microRNAs. Rna 2012, 18(3):472-484. 

76. Piriyapongsa J, Marino-Ramirez L, Jordan IK: Origin and evolution of human 
microRNAs from transposable elements. Genetics 2007, 176(2):1 323-1 337. 

77. Smalheiser NR, Torvik VI: Mammalian microRNAs derived from genomic 
repeats. Trends Genet 2005, 21 (6):322-326. 

78. Piriyapongsa J, Jordan IK: A family of human microRNA genes from 
miniature inverted-repeat transposable elements. PLoS One 2007, 2(2):e203. 

79. Lehnert S, Van Loo P, Thilakarathne PJ, Marynen P, Verbeke G, Schuit FC: 
Evidence for co-evolution between human microRNAs and Alu-repeats. 
PLoS One 2009, 4(2):e4456. 

80. Ono M, Scott MS, Yamada K, Avolio F, Barton GJ, Lamond Al: Identification 
of human miRNA precursors that resemble box C/D snoRNAs. Nucleic 
Acids Res 201 1, 39(9):3879-3891. 

81. Scott MS, Avolio F, Ono M, Lamond Al, Barton GJ: Human miRNA 
precursors with box H/ACA snoRNA features. PLoS Connput Biol 2009, 
5(9):e1 000507. 



82. Brameier M, Herwig A, Reinhardt R, Walter L, Gruber J: Human box C/D 
snoRNAs with miRNA like functions: expanding the range of regulatory 
RNAs. Nucleic Acids Res 201 1, 39(2):675-686. 

83. Grundhoff A, Sullivan CS: Virus-encoded microRNAs. Virology 201 1, 
411(2):325-343. 

84. Singh J, Singh CP, Bhavani A, Nagaraju J: Discovering microRNAs from 
Bombyx mori nucleopolyhedrosis virus. Virology 2010, 407(1 ):1 20-1 28. 

85. Singh CP, Singh J, Nagaraju J: A baculovirus-encoded MicroRNA (miRNA) 
suppresses its host miRNA biogenesis by regulating the exportin-5 
cofactor Ran. J Virol 2012, 86(15):7867-7879. 

86. Guan DG, Liao JY, Qu ZH, Zhang Y, Qu LH: mirExplorer: detecting 
microRNAs from genome and next generation sequencing data using 
the AdaBoost method with transition probability matrix and combined 
features. RNA Biol 201 1, 8(5):922-934. 

87. Zywicki M, Bakowska-Zywicka K, Polacek N: Revealing stable processing 
products from ribosome-associated small RNAs by deep-sequencing 
data analysis. Nucleic Acids Res 2012, 40(9):40 13-4024. 

88. Taft RJ, Glazov EA, Lassmann T, Hayashizaki Y, Carninci P, Mattick JS: Small 
RNAs derived from snoRNAs. Rna 2009, 15(7):1 233-1 240. 

89. Ender C, Krek A, Friedlander MR, Beitzinger M, Weinmann L, Chen W, Pfeffer 
S, Rajewsky N, Meister G: A human snoRNA with microRNA-like functions. 
Mol Cell 2008, 32(4):5 19-528. 

90. Ma Z, Xue Z, Zhang H, Li Y, Wang Y: Local and global effects of Mg2-H on 
Ago and miRNA-target interactions. J Mol Model 2012, 18(8):3769-3781. 

91 . Cao J, Tong C, Wu X, Lv J, Yang Z, Jin Y: Identification of conserved 
microRNAs in Bombyx mori (silkworm) and regulation of fibroin L chain 
production by microRNAs in heterologous system. Insect Biochenn Mol Biol 

2008, 38(1 2):1 066-1 071. 

92. Saxena SK, Rybak SM, Davey RT Jr, Youle RJ, Ackerman EJ: Angiogenin is a 
cytotoxic, tRNA-specific ribonuclease in the RNase A superfamily. J Biol 
Chem 1992, 267(30):21 982-21 986. 

93. Yamasaki S, Ivanov P, Hu GF, Anderson P: Angiogenin cleaves tRNA and 
promotes stress-induced translational repression. J Cell Biol 2009, 
185(1):35-42. 

94. Fu H, Feng J, Liu Q, Sun F, Tie Y, Zhu J, Xing R, Sun Z, Zheng X: Stress 
induces tRNA cleavage by angiogenin in mammalian cells. FEBS Lett 

2009, 583(2):437-442. 

95. Chung WJ, Okamura K, Martin R, Lai EC: Endogenous RNA interference 
provides a somatic defense against Drosophila transposons. Curr Biol 
2008, 18(11):795-802. 

96. Ghildiyal M, Seitz H, Horwich MD, Li C, Du T, Lee S, Xu J, Kittler EL, Zapp ML, 
Weng Z, et al: Endogenous siRNAs derived from transposons and mRNAs 
in Drosophila somatic cells. Science 2008, 320(5879):1 077-1 081. 

97. Watanabe T, Totoki Y, Toyoda A, Kaneda M, Kuramochi-Miyagawa S, Obata 
Y, Chiba H, Kohara Y, Kono T, Nakano T, et al: Endogenous siRNAs from 
naturally formed dsRNAs regulate transcripts in mouse oocytes. Nature 
2008, 453(71 94):539-543. 

98. Nie Z, Xu J, Chen J, Lv Z, Wang D, Sheng Q, Wu Y, Wang X, Wu X, Zhang Y: 
Expression analysis and characteristics of profilin gene from silkworm, 
Bombyx mori. AppI Biochem Biotechnol 2009, 158(1):59-71. 

99. Gomi S, Majima K, Maeda S: Sequence analysis of the genome of Bombyx 
mori nucleopolyhedrovirus. J Gen Virol 1999, 80(Pt 5):1 323-1 337. 

100. Li R, Yu C, Li Y, Lam TW, Yiu SM, Kristiansen K, Wang J: S0AP2: an 
improved ultrafast tool for short read alignment. Bioinformatics 2009, 
25(1 5):1 966-1 967. 

101. Langmead B, Trapnell C, Pop M, Salzberg SL: Ultrafast and memory- 
efficient alignment of short DNA sequences to the human genome. 
Genome Biol 2009, 10(3):R25. 

102. Huang Y, Niu B, Gao Y, Fu L, Li W: CD-HIT Suite: a web server for 
clustering and comparing biological sequences. Bioinformatics 2010, 
26(5):680-682. 

103. Yang L, Lu X, Liu Y, Lv Z, Chen J, Yu W, Zhang Y, Nie Z: Expression analysis 
of miRNAs in BmN cells. Gene 2012, 505(2):240-245. 



doi:1 0.1 186/1471-2164-14-661 

Cite this article as: Nie et al.: RIP-seq of BmAgo2-associated small RNAs 
reveal various types of small non-coding RNAs in the silkworm, Bombyx 
mori. BMC Genomics 2013 14:661. 



